Method and apparatus for relative navigation using reflected gps signals

ABSTRACT

A method and system to passively navigate an orbiting moving body towards an orbiting target using reflected GPS signals. A pair of antennas is employed to receive both direct signals from a plurality of GPS satellites and a second antenna to receive GPS signals reflected off an orbiting target. The direct and reflected signals are processed and compared to determine the relative distance and position of the orbiting moving body relative to the orbiting target.

GOVERNMENT INTEREST

The embodiments of the invention described herein were made by employeesof the United States Government, and may be manufactured and used by orfor the United States Government for governmental purposes withoutpayment of any royalties thereon or therefor.

BACKGROUND OF THE INVENTION

1. Field of the Invention

The present invention relates to navigating a moving body with the useof GPS signals and more particularly to the use of reflected GPS signalsto passively navigate an orbiting body towards an orbiting target.

2. Description of the Related Art

Spacecraft require equipment that calculates the orbital position,velocity, time, and attitude. The advent of the Global PositioningSystem (“GPS”), increased the desirability of using onboard processingsystems for spacecraft position determination.

The GPS satellite system is a collection of satellites that can be usedfor missile, satellite, aircraft, and terrestrial navigation. Each GPSsatellite broadcasts its own ephemeris and time, thereby allowing a GPSreceiver to determine its position. Typically, a GPS receiver calculatesits position from the simultaneous observation of any four GPSsatellites in order to calculate its position. A civil grade GPSreceiver can accurately determine its position within 20 meters,determine its velocity within 0.6 meters/second, and the current timewith a 10-nanosecond accuracy. At least twenty-four valid GPS satellitesare in the GPS constellation at any given time.

A GPS receiver is an autonomous instrument that transforms signals fromGPS satellites into point solutions for spacecraft navigation. CurrentGPS receivers have a radio frequency section for receiving andconverting signals received from a spacecraft's antennas. The digitizedsignals are then forwarded to one or more correlators, controlled by thereceivers own processor. The correlators look for matches between theincoming signal and the C/A code for different satellites. When asatellite lock occurs, or the incoming signal matches an internallygenerated pseudo random noise (PRN) code, the receiver's processor isnotified. The processor contains executable code to generate apseudo-range or a line-of-sight distance to the satellite. The processormay also contain the executable code of an orbit propagator. Thepseudo-range is a measurement input to the navigation filter thatcalculates point solutions for determining the orbit of the spacecraft.While such GPS systems may be employed to position and track an orbitingmoving body such as the Space Shuttle, such systems do not provide theability to determine the position relative to another orbiting body.Conventional GPS receivers can not passively provide relative navigationor relative position data which are necessary during formation flying,docking, or other approach of space craft to proximity to anotherorbiting body.

The space environment poses additional obstacles to that of land basedtracking systems. Because of the speeds associated with spacecraftmotion, GPS satellites are in view for much shorter time periods thanwhen viewed from a ground reference. The appearance and disappearance ofGPS satellites from view causes their output signals to slew through theentire range of the 45-kilohertz Doppler shift. Typical space craftformation flying and autonomous rendezvous missions rely on activetransmissions schemes such as RADAR (RAdio Detection and Ranging) andLIDAR (Light Detection and Ranging) which require specialized hardwarerequiring additional mass and power consumption. With the extreme weightand power restraints associated with space flight, there is a particularneed for a reliable navigation/positioning system which reduces weightand power consumption requirements, as well as taking advantage of freeGPS signals.

The use of Reflected GPS signals per se is also known. U.S. patentapplication publication 2006/0077094 entitled Information “InformationGathering Using reflected Signals” discloses a method of using reflectedsignals in a land based application to monitor the presence or absenceof an object such as if a vehicle occupies a parking space and is herebyincorporated herein in its entirely. This system is not only ill suitedfor orbital applications, but does not make use of the reflected signalto navigate to a target.

It is also known to use Bistatic Radar Systems where a radar transmitterlocated remotely from the projectile (such as onboard a ship)illuminates the target and the reflected returns are received by areceiver located on the projectile. The tracking data from the radarmeasurements are then used to calculate the proper guidance signals todirect the projectile to the target. An example of a bistatic radarsystem is discloses ill U.S. Pat. No. 6,653,972 which in incorporatedherein it is entirely. These systems require the use and control of anactive radar transmitter and do not provide a suitable globalpositioning function. Rather, these systems are employed to simplyacquire and hone in on a particular target. Therefore, such Bistaticradar systems are not suited for orbital positioning and navigation orhave the ability to rely on a passive autonomous system.

One of the objects of the present invention is to utilize reflected GPSsignals to provide passive autonomous relative navigation of an orbitingspacecraft towards an orbiting body.

SUMMARY OF THE INVENTION

The present invention is directed to a method and system to passivelynavigate an orbiting moving body towards an orbiting target usingreflected GPS signals. A pair of receive antennae are employed. One toreceive direct signals from a plurality of GPS satellites and a secondantenna to receive GPS signals reflected off an orbiting target. Thedirect and reflected signals are processed and compared to determine therelative distance and position of the orbiting moving body relative tothe orbiting target.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is a diagram depicting the relative geometry between an orbitingspacecraft, target and satellite.

FIGS. 2-3 are diagrams depicting the earth centered inertial coordinateframe.

FIG. 4 is a diagram depicting the visibility of GPS signals to a body inorbit.

FIG. 5 is a block diagram of a system to receive and process reflectedsignals according to a preferred embodiment of the present invention.

FIG. 6 is a flow diagram depicting the algorithm to receive and processreflected signals according to a preferred embodiment of the presentinvention.

DETAILED DESCRIPTION OF A PREFERRED EMBODIMENT

A preferred embodiment of the invention and the various features andadvantageous details thereof are explained more fully with reference tothe accompanying drawings and detailed in the following description. Itshould be noted that the features illustrated in the drawings are notnecessarily drawn to scale and descriptions of well-known components andprocessing techniques are omitted so as to not unnecessarily obscure theembodiments of the invention. The examples used herein are intendedmerely to facilitate an understanding of ways in which the embodimentsof the invention may be practiced and to further enable those of skillin the art to practice the embodiments of the invention. Accordingly,the examples should not be construed as limiting the scope of theembodiments of the invention.

The present invention is directed to a method and apparatus forpassively navigating an orbiting space craft towards and orbiting targetusing reflected GPS signal. FIG. 1 depicts the relative geometry betweenthe space shuttle, Hubble space telescope and a GPS satellite. Bycomparing both signals received directed from the GPS satellite andsignals reflected off the Hubble, a navigation system employed on thespace shuttle can passively and autonomous determine its positionrelative to the Hubble thereby providing navigational aid passively andautonomously. Acquisition of and tracking of GPS signals are well knownin the art and need not be discussed in detail. U.S. published patentapplication US 2006-0082496 to Winternitz et al. discloses aradiation-hardened fast acquisition/weak signal tracking system andmethod and is hereby incorporated herein by referenced in its entirety.

The actual GPS measurements recorded from any of the constellation ofsatellites are the pseudorange and psuedorange rate. The pseudorange, ρ,is based on the time of flight between the GPS receiver and the GPStransmit antenna. This is used to form a range measurement, and thussolve for a position fix when enough satellites are available. The othermeasurement available is pseudorange rate that is essentially the rateof change along the line of sight vector to the satellite. Thepseudorange rate, dotρ, can also be viewed as the Doppler shift on thesignal. The pseudorange rate equation is equivalent to the carrier rate,Φ which tends to be a less noisy measurement.

The measured pseudorange is the difference between the system time atwhich the signal left the GPS satellite and the system time at which thesignal reached the user. This is given by

ρ=r+c(t _(u) +δt+δt _(D))

r is the geometric range is given by r=c(T_(u)-T_(s)), which is theexact time of flight between the user and the GPS satellite. t_(u) isthe receiver offset from system time, δ_(t) is the offset of thesatellite from system time, and δ_(tD) is due to other measurementerrors.

The GPS constellation is segmented into 3 parts, the user, the control,and space segment. The three main sources for the space segment aresatellite clock stability, clock perturbations, and selectiveavailability. For the control segment, it's primarily ephemerisprediction error. For the user segment it is ionospheric delay,multipath, receiver noise, and resolution.

The GPS satellite clock error is the residual error from a second orderpolynomial fit, adjusting for relativistic effects, the clock bias,clock drift, and frequency drift. The Satellite Clock error is on theorder of 3.0 m. The satellite ephemeris error is caused from thedeviation in the current satellite ephemeride stored onboard the GPSsatellite with the satellite's true position. These ephemeris values areused to form the line of sight vectors from each SV to the user, thuseach perturbation adds error into the pseudorange. This is normally onthe order of 4.2 m for a typical ground based GPS receiver.

The GPS system also experiences relativistic effects due to thegravitational potential differences. Space users also experienceionospheric effects. This is generally seen as a bias in themeasurement. The ionospheric effect is modeled off of the path length ofthe signal traveled through the ionosphere. There are also measurementerrors from within the receiver. These are from the tracking loopsinternal to the GPS receiver. The dominant sources of error introducedhere are the thermal noise jitter and dynamic stress error. Thesecondary sources are given by the hardware and software resolution, andclock drift. The last sources are multipath and shadowing effects. Forthe purposes of rendezvous, we are actually interested in tracking thesemultipath sources in order to solve for a relative position. The GPSpseudorange error budget is on the order of 8.0 m for a standard GPSreceiver located on the ground.

The direct pseudorange rate is a measurement of the rate of change ofthe psuedorange between the user and the GPS satellite. It is given by

{dot over (ρ)}=(V _(u) −V _(a))·lôs+ε

That is the rate of change of the psuedorange along the line of sightvector and the user. The error is mainly determined by the clock driftcomponent. Generally these measurements are noisy compared to that ofthe pseudorange.

The constructed measurements for the bistatic ranging problem are thedifferenced reflected psuedorange and the reflected psuedorange rate.The reflected psuedorange, ρ_(R) is differenced with the direct signalpsuedorange to remove the common mode errors from the measurement. Someof the common mode errors removed are front-end jitter, receiver bias,thermal noise, and line bias, leaving us primarily with ionosphericeffects, multipath and jitter error contributions. The psuedorange rateis the additive Doppler shifts along the line of sight vectors betweenHubble and the Shuttle, and between Hubble and the GPS constellation.

This difference reflected pseudorange measurement is constructed fromthe difference in the direct path signal and the reflected signal. Theequation is given by:

ρ_(ρr)=ρ_(R)−ρ_(D)

Expressing as a geometric relation we arrive at:

ρ_(ρr) =r _(sh) +r _(hg) −r _(gs)+ε

That is the sum of the geometric range from the Shuttle to Hubble,Hubble to the GPS satellite, and the Shuttle to the GPS satellite, plusthe residual errors.

The psuedorange rate can be seen as a Doppler shift on the signal. Thereflection can be viewed as a retransmission, making it an additiveprocess.

{dot over (ρ)}_(r)=(V _(s) −V _(h))lôs_(sh)+(V _(h) −V _(g))lôs_(hg)

Where the lôs_(sh) is the line of sight vector between the Shuttle andHubble, and lôs_(sh) is the line of sight between Hubble and a GPSsatellite. The primary sources of error contributing to the reflectedpsuedorange rate are the same of that as the direct measurement,however, the quality of the incoming signal is degraded due to thereflective process and the actual effect is unknown.

Typical GPS receivers use a simple least squares algorithm to determineposition and time. It requires at least 4 visible satellites to solvefor the four unknowns x,y,z, and the user clock offset from the GPSconstellation. In addition to position estimates, some also do perform afinite difference of the positions to solve for velocity as well.

In more advanced receivers, a Kalman filter may be implemented in orderto utilize the dynamic information, as well as the measurements toimprove the accuracy of the receiver. On orbit, the GPS receiver willexperience large dynamic stress and motion with respect to the GPSsatellites, it is almost necessary to filter this heavily to provide agood navigation solution. The method of the present invention implementsan extended Kalman filter which will now be described.

The Kalman filter represents an optimal filter for a linear system withadditive gaussian noise with zero mean. Since the dynamics of anorbiting body are well understood and their equations are wellformulated, it is possible to exploit that information when forming theKalman filter when applied to the relative navigation problem. Thecanonical form of the Kalman filter is the optimal filter for a linearsystem, incorporating a blending of known dynamics, with that of themeasurements. The discrete Kalman filter can be broken up into twodistinct steps, the time update and measurement update.

The dynamical equations in a standard linear system are in the form:

{dot over (X)}=Ax+Bu+d

y=Hx+n

These equations are continuous and need to be discretized in order toimplement them in the Kalman filter.This discretization is given by:

x _(k+1)=Φ_(k) x+Γ _(k) u _(k) +d _(k)

y _(k) =Hx _(k) +n _(k)

This relates timestep k−1, to timestep k. In order to calculate thestate evolution between timesteps, for a linear system, the statetransition matrix Φ is given by:

Φ_(k)=e^(AΔT)

Where

ΔT=t _(k+1) −t _(k).

The evolution of Γ is a bit more complex, it is given by:

Γ_(k) = ∫₀^(T)^(A Δτ)B u(τ)τ

This model only concerns the unforced dynamics, so without loss ofgenerality, B is assigned 0 and Γ is dropped. Now these equations can beapplied to the Kalman filter form.

The first step in the Kalman filter is to propagate the measurement fromtime-step k to the time-step k+1. This is done by multiplying by thestate transition matrix, or integrating from t_(k) to t_(k+1).

{circumflex over (x)} _(k+1)=Φ_(k) {circumflex over (x)} _(k)+Γ_(k) u_(k)

Propagating the error covariance matrix.P is defined as

P=E[({circumflex over (x)}−x)({circumflex over (x)}−x)^(T)]

it represents the variance in the estimator errors.The error covariance propagation step is given as follows:

P _(k+1)=Φ_(k) P _(k)Φ_(k) ^(T) +Q

The first step of the measurement update is to form the Kalman Gain,K_(K). The Kalman gain is the optimal estimate to minimize the errorsalong the diagonal of the covariance matrix.

K _(k) =P _(k) H _(k) ^(T)(H _(k) P _(k) H _(k) ^(T) +R)⁻¹

The next step is to apply the gain to the measurement error. The Hmatrix is the measurement connection matrix, it connects the states withthe measurement at time-step k.

{circumflex over (x)} _(k) ={circumflex over (x)} _(k) +K _(k)(z _(k) −H_(k) {circumflex over (x)} _(k) )

The final step is to correct the error covariance for time-step k.

P _(k)=(I−K _(k) H)P _(k)

While die equations are relatively straight forward in the discreteKalman filter, they require a measurement model and the plant dynamicsare linear systems. Unfortunately, in the real world most systems arenonlinear. In particular, GPS has a non-linear measurement model.

The solution to a nonlinear measurement model can be done byapproximating Y_(K). Given a nonlinear measurement model Y=h(x)

${{{y_{k} = {H_{k}x_{k}}}{H_{k} = \frac{\partial h}{\partial x}}}}_{{\hat{x}}_{k}}$

This approximation is called an extended Kalman filter.The last two parts of a Kalman filter are the measurement covariance,and process noise matrices, R and Q respectively. The measurementcovariance matrix consists of the variance of noise on the signal. Fromthe canonical form:

y _(k) =H _(k) x _(k) +n _(k)

it can be shown that

R=E[n_(k)n_(k) ^(T)]

Assuming uncorrelated gaussian white noise with zero mean, R is givenby:

$R = \begin{bmatrix}{\sigma_{\rho}^{2}I_{nxn}} & 0_{nxn} \\0_{nxn} & {\sigma_{\overset{.}{\rho}}^{2}I_{nxn}}\end{bmatrix}$

The process noise matrix is given by the variance un-modeled forces inthe dynamic model. Assuming we have a model:

x _(k)=Φ_(k) x _(k)+Γ_(k) u _(k) +d _(k)

noise actually enters the system as a disturbance on the input.

Q=Γ_(k)E[d_(k)d_(k) ^(T)]Γ_(k) ^(T)

Unfortunately, the equation to describe the matrix Q is difficult toobtain because it requires calculation of Γ. In the simulations, Q istreated as a free variable, adjusting it to give the optimalperformance, and reduce the process covariance bounds such that theywere roughly equal to the measured variance in R.

In the case of relative navigation, it is an essential to be able tosolve for your relative position as well as possible. One of the bestways is to improve your solution is to have an accurate dynamic model.This way, one can rely on the dynamic model to propagate throughmeasurement outages and noisy measurements. This is especially true inthe case of applying reflected GPS measurements.

Defining the necessary components to construct the filter is necessarybefore beginning formulation of the Kalman filters for the relativenavigation problem. Given a dynamical system:

{dot over (x)}=f(x)

y=h(x)

It is thus necessary to define a state vector, X. The state vectors usedin an absolute model, are the absolute positions and velocities of theShuttle and HST, as well as the bias state for the navigation solution.Typically a GPS receiver would include the bias drift as anothercomponent. For this simulation, since the rendezvous is over a shortperiod of time, bias drift can be neglected.

x=[r_(s),v_(s),r_(h),v_(h),b]^(T)

For a completely relative model, Hill's equations are used to propagatethe states. This requires formulating the problem entirely in a relativeframe. The states in this case are:

X=[Δ_(r),Δ_(v)]^(T)

In this formulation, Δ_(r) is the relative position in hill's frame andΔ_(v) is the velocity. Typically the dynamic model is a nonlinearequation given by a f(dot)=f(x). In this analysis two different modelare examined. The first model is based purely on the 2 body equation ofmotion given by Hill's equations. This model is a linearizedrepresentation of the relative spacecraft dynamics and has the form

{dot over (x)}=Ax.

where A is given by the partials of f(x) evaluated along a referencetrajectory.

${{A = \frac{\partial f}{\partial x}}}_{{\hat{x}}_{k}}$

Hill's model is conveniently linear, simplifying the application of theKalman filter. Hill's model is only valid for close proximity maneuversand requires a special coordinate frame to operate in, centered on acircular orbit for one of the vehicles.

Hill's equations of motion are derived in the {circumflex over (R)}, Ŝ,Ŵ coordinate frame as depicted in FIG. 2.

{circumflex over (R)} is the radial direction from the center of theearth,Ŝ is the projection of the velocity in the plane perpendicular to{circumflex over (R)}, and Ŵ is the direction of angular momentumperpendicular to the orbital plane.The transformations from the IJK coordinate frame are as follows:

$\hat{R} = \frac{R}{R}$ Ŝ = Ŵ × R̂$\hat{W} = \frac{\hat{R} \times \overset{\overset{.}{\hat{}}}{R}}{{\hat{R} \times \overset{\overset{.}{\hat{}}}{R}}}$

Where r and {dot over (r)} are the position and velocity vectors in theIJK reference frame. Since Hill's frame is a relative frame, typicallyone spacecraft is chosen to be the reference orbit. This orbit isideally circular. In this case, the orbit chosen is the space shuttle.This is denoted by r_(tgt), the other spacecraft is typically called theinterceptor or r_(int). Hill's equations, also known as the ClohessyWiltshire equations which are the governing forces for relative two-bodymotion. These equations are well known and have an explicit solution.

{umlaut over (x)}−2w{dot over (y)}−3w ² x+fx=0

ÿ+2w{dot over (x)}+fy=0

{umlaut over (z)}+w ² z+fz=0

or, by definingx^(T)=(x,{dot over (x)}), then provides the form {dot over (x)}=f(x).Now, taking partials with respect to the states, arriving at thecanonical linear system form {dot over (x)}=Ax.

$\overset{.}{x} = {{\begin{bmatrix}0 & 0 & 0 & 1 & 0 & 0 \\0 & 0 & 0 & 0 & 1 & 0 \\0 & 0 & 0 & 0 & 0 & 1 \\{3\; \omega^{2}} & 0 & 0 & 0 & {2\omega} & 0 \\0 & 0 & 0 & {{- 2}\omega} & 0 & 0 \\0 & 0 & \omega^{2} & 0 & 0 & {0\;}\end{bmatrix}x} + {\begin{bmatrix}0 & 0 & 0 \\0 & 0 & 0 \\0 & 0 & 0 \\1 & 0 & 0 \\0 & 1 & 0 \\0 & 0 & 1\end{bmatrix}u}}$

Hill's equations being in a linear form fit well into the linear systemstoolset. If we assume a constant A over each discrete time step we canapproximate:

Φ(t_(k + 1,)t_(k)) = ^((A Δ t))$A = \left. \frac{\partial f}{\partial x} \right|_{{\hat{x}}_{k}}$

This formulation provides a state transition matrix in linear form thatcan be easily applied to both the state propagation steps and thecovariance propagation. Although the relative dynamics are minimal incomparison between spacecraft, their actual perturbations on thedynamics are quite significant in low earth orbit. In addition toincreasing the fidelity of the dynamic model, it allows customization ofthe orbits and adds other effects such as thruster firing, anddifferences in their drag coefficients.

The coordinate frame used in the high fidelity dynamic model is theEarth Centered Inertial frame. This reference frame was chosen becauseit greatly simplifies the equations of motion. The ECI coordinate systemis also known as the IJK system. Referring to FIG. 3 the Î axis ispointed towards the vernal equinox, the Ĵ axis is 90° to the east alongthe equatorial plane, and the {circumflex over (K)} axis extends throughthe North Pole.

The full nonlinear dynamic equations are given by:

${\overset{¨}{r}}_{x} = {{- {\frac{\mu \; x}{r^{3}}\left\lbrack {1 - {\frac{3}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{\frac{5}{2}\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} - {\frac{1}{2}C_{D}\frac{A}{m}\rho \; {V\left( {\overset{.}{x} + {\omega_{E}y}} \right)}}}$${\overset{¨}{r}}_{y} = {{- {\frac{\mu \; y}{r^{3}}\left\lbrack {1 - {\frac{3}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{\frac{5}{2}\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} - {\frac{1}{2}C_{D}\frac{A}{m}\rho \; {V\left( {\overset{.}{y} + {\omega_{E}x}} \right)}}}$${\overset{¨}{r}}_{z} = {{- {\frac{\mu \; z}{r^{3}}\left\lbrack {1 - {\frac{3}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{\frac{5}{2}\frac{z^{2}}{r^{2}}} - 3} \right)}} \right\rbrack}} - {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\overset{.}{\; z}}}$

The J₂ term is the zonal harmonic constant given by 1082.7×10⁻⁶. ω_(E)is the rotation rate of the earth, R_(E) is the radius of the earth, mis the objects mass, A is the cross sectional area, and V is therelative wind. The relative wind is given by the equation

$V = \sqrt{{\overset{.}{r}}_{x}^{2} + {2\omega_{E}{\overset{.}{r}}_{x}r_{y}} + {\overset{.}{r}}_{y}^{2} - {2\omega_{E}{\overset{.}{r}}_{y}r_{x}} + {\omega_{E}^{2}r_{y}^{2}} + {\overset{.}{r}}_{z}^{2}}$

The value of ρ is given by an exponential model derived in:

ρ=ρ₀ e ^((−α−α) ^(o) ^()/H)

Where ρ₀=3.275e−12 kg/m³, α is the current altitude, α₀ is the referencealtitude, 400000 m and H is the scaling height 58515 m. Using the samefrom the simple GPS Kalman filter, choosing r as position and v asvelocity, the equations it is now appropriate to write these equationsin terms of first order derivatives, to fit the canonical form, are

{dot over (r)}=v and {dot over (v)}=g(r,v)

In this formulation, g is the nonlinear force model for the orbitaldynamics. Now, let the state vector x=┌r,v┐^(T).

{dot over (x)}=[v,g(r,v)]^(T)

or

${f(x)} = \begin{bmatrix}v_{x} \\v_{y} \\v_{z} \\{\overset{¨}{r}}_{x} \\{\overset{¨}{r}}_{y} \\{\overset{¨}{r}}_{z}\end{bmatrix}$

The equations are in the appropriate form to apply the Kalman Filter.Since the Kalman filter's covariance propagation is a linear process, itis necessary to formulate the partial derivatives with respect to thestate vector. Although the partials are not used to propagate thedynamics, they are needed when forming the state transition matrix forthe covariance updates.

$\mspace{79mu} {\frac{\partial{\overset{.}{r}}_{x}}{\partial\overset{.}{x}} = {{1\mspace{14mu} \frac{\partial{\overset{.}{r}}_{y}}{\partial\overset{.}{y}}} = {{1\mspace{14mu} \frac{\partial{\overset{.}{r}}_{z}}{\partial\overset{.}{z}}} = 1}}}$$\frac{\partial{\overset{¨}{r}}_{x}}{\partial x} = {{- {\frac{\mu}{r^{3}}\left\lbrack {1 - {\frac{3}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{\frac{5}{2}\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} + {3{\frac{\mu \; x^{2}}{r^{5}}\left\lbrack {1 - {\frac{5}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{7\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{x\left( {\overset{.}{x} + {\omega_{E}y}} \right)}{r\; H}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{{x\left( {\overset{.}{x} + {\omega_{E}y}} \right)}\left( {\overset{.}{y} - {\omega_{E}x}} \right)}{V}}}$$\frac{\partial{\overset{¨}{r}}_{x}}{\partial y} = {{3{\frac{\mu \; x\; y}{r^{5}}\left\lbrack {1 - {\frac{5}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{7\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{y\left( {\overset{.}{x} + {\omega_{E}y}} \right)}{r\; H}} - {\frac{1}{2}C_{D}\frac{A}{m}{\rho \;\left\lbrack {{V\; \omega_{E}} + \frac{{\omega_{E}\left( {\overset{.}{x} + {\omega_{E}y}} \right)}^{2}}{V}} \right\rbrack}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{x}}{\partial z} = {{3{\frac{\mu \; x\; z}{r^{5}}\left\lbrack {1 - {\frac{5}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{7\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{z\left( {\overset{.}{x} + {\omega_{E}y}} \right)}{r\; H}}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{x}}{\partial\overset{.}{x}} = {{- \frac{1}{2}}C_{D}\frac{A}{m}{\rho \;\left\lbrack {V + \frac{\left( {\overset{.}{x} + {\omega_{E}y}} \right)^{2}}{V}} \right\rbrack}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{x}}{\partial\overset{.}{y}} = {{- \frac{1}{2}}C_{D}\frac{A}{m}\rho \frac{\left( {\overset{.}{y} - {\omega_{E}x}} \right)\left( {\overset{.}{x} - {\omega_{E}y}} \right)}{V}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{x}}{\partial\overset{.}{z}} = {{- \frac{1}{2}}C_{D}\frac{A}{m}\rho \; \frac{\overset{.}{z}\left( {\overset{.}{x} + {\omega_{E}y}} \right)}{V}}}$$\frac{\partial{\overset{¨}{r}}_{y}}{\partial x} = {{3{\frac{\mu \; x\; y}{r^{5}}\left\lbrack {1 - {\frac{5}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{7\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{x\left( {\overset{.}{y} + {\omega_{E}x}} \right)}{r\; H}} - {\frac{1}{2}C_{D}\frac{A}{m}{\rho \;\left\lbrack {{V\; \omega_{E}} + \frac{{\omega_{E}\left( {\overset{.}{y} + {\omega_{E}y}} \right)}^{2}}{V}} \right\rbrack}}}$$\frac{\partial{\overset{¨}{r}}_{y}}{\partial y} = {{- {\frac{\mu \;}{r^{3}}\left\lbrack {1 - {\frac{3}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{5\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} + {3{\frac{\mu \; y^{2}}{r^{5}}\left\lbrack {1 - {\frac{5}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{7\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{y\left( {\overset{.}{y} + {\omega_{E}x}} \right)}{r\; H}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; \frac{{\omega_{E}\left( {\overset{.}{y} + {\omega_{E}x}} \right)}\left( {\overset{.}{x} - {\omega_{E}y}} \right.}{V}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{y}}{\partial z} = {{3{\frac{\mu \; x\; z}{r^{5}}\left\lbrack {1 - {\frac{5}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{7\frac{z^{2}}{r^{2}}} - 1} \right)}} \right\rbrack}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{z\left( {\overset{.}{y} + {\omega_{E}y}} \right)}{r\; H}}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{y}}{\partial\overset{.}{x}} = {{- \frac{1}{2}}C_{D}\frac{A}{m}\rho \; \frac{\left( {\overset{.}{y} - {\omega_{E}x}} \right)\left( {\overset{.}{x} - {\omega_{E}y}} \right)}{V}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{y}}{\partial\overset{.}{y}} = {{- \frac{1}{2}}C_{D}\frac{A}{m}{\rho \left\lbrack {V + \; \frac{\left( {\overset{.}{y} + {\omega_{E}x}} \right)^{2}}{V}} \right\rbrack}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{y}}{\partial\overset{.}{z}} = {{- \frac{1}{2}}C_{D}\frac{A}{m}\rho \; \frac{\overset{.}{z}\left( {\overset{.}{y} + {\omega_{E}z}} \right)}{V}}}$$\frac{\partial{\overset{¨}{r}}_{z}}{\partial x} = {{3{\frac{\mu \; x\; z}{r^{5}}\left\lbrack {1 - {\frac{5}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{7\frac{z^{2}}{r^{2}}} - 3} \right)}} \right\rbrack}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{\overset{.}{z}x}{r\; H}} - {\frac{1}{2}C_{D}\frac{A}{m}\rho \; \frac{\omega_{E}{\overset{.}{z}\left( {\overset{.}{y} + {\omega_{E}x}} \right)}}{V}}}$$\frac{\partial{\overset{¨}{r}}_{z}}{\partial z} = {{3{\frac{{\mu y}\; z}{r^{5}}\left\lbrack {1 - {\frac{5}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{7\frac{z^{2}}{r^{2}}} - 3} \right)}} \right\rbrack}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{\overset{.}{z}y}{r\; H}} - {\frac{1}{2}C_{D}\frac{A}{m}{\rho\left( \; \frac{\omega_{E}{\overset{.}{z}\left( {\overset{.}{x} - {\omega_{E}y}} \right)}}{V} \right)}}}$$\frac{\partial{\overset{¨}{r}}_{z}}{\partial z} = {{- {\frac{\mu \;}{r^{3}}\left\lbrack {1 - {\frac{3}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{5\frac{z^{2}}{r^{2}}} - 3} \right)}} \right\rbrack}} + {3{\frac{\mu \; y^{2}}{r^{5}}\left\lbrack {1 - {\frac{5}{2}{J_{2}\left( \frac{R_{E}}{r} \right)}^{2}\left( {{7\frac{z^{2}}{r^{2}}} - 5} \right)}} \right\rbrack}} + {\frac{1}{2}C_{D}\frac{A}{m}\rho \; V\frac{z\overset{.}{z}}{r\; H}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{z}}{\partial\overset{.}{x}} = {{- \frac{1}{2}}C_{D}\frac{A}{m}\rho \; \frac{\overset{.}{z}\left( {\overset{.}{x} + {\omega_{E}y}} \right)}{V}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{z}}{\partial\overset{.}{y}} = {{- \frac{1}{2}}C_{D}\frac{A}{m}\rho \; \frac{\overset{.}{z}\left( {\overset{.}{y} + {\omega_{E}x}} \right)}{V}}}$$\mspace{79mu} {\frac{\partial{\overset{¨}{r}}_{z}}{\partial\overset{.}{z}} = {{- \frac{1}{2}}C_{D}\frac{A}{m}{\rho \left( {V + \; \frac{{\overset{.}{x}}^{2}}{V}} \right)}}}$

Hill's equations of motions are conveniently in linear form:x_((dot))=A_(x)+B_(u). The state transition matrix can simply beapproximated by the matrix exponential

e^(AΔT).

Where A is the linear system A matrix, and T is the time interval. Inthis particular implementation the time step is on the order of 1second, thus T=1 This is both used in the dynamics, and the covariancepropagation. In the full nonlinear case estimates would diverge after afew iterations using a first order approximation of the system. Thus thestate is integrated over time steps of one second. The covariancematrices are propagated linearly, using,

P _(k+1)=Φ(t _(k+1) ,t _(k))P _(k)Φ(t _(k+1) ,t _(k))^(T) +Q.

Solving for Φ in a given linear system, we have:

{dot over (Φ)}(t _(k+1) ,t _(k))=A(t _(k))Φ(t _(k))

Φ(t _(k+1) ,t _(k))=∫_(t) _(k) ^(t) ^(k+1) A(τ)Φdτ

and arrive at an approximate

Φ(t_(k−1),t_(k)).

Since A is not constant over the time interval t_(k+1), t_(k), we cannotrely on the matrix exponential approximation for this integration.In order to incorporate the variation of the states between time steps,the added dynamic equation {dot over (Φ)}=AΦ is included within thedynamic model propagation. This allows integrating the dynamics andsolving for φ between time steps, while varying A with the current stateestimates. In order to solve for φ, the initial condition vector of thestates and the identity reshaped into a vector is required. Theresultant φ at the end of the integration is then used in the linearpropagation of P.

Implementing the filter with Hill's model is fairly straightforward.Since Hill's model is in the relative frame, the shuttle's position isimplicitly included in the states, any parameters derived from theposition, can be assumed to be provided by the GPS receiver's navigationsolution.

The equation describing the reflected Pseudorange Measurement is givenby:

h(X)=ρ_(r) =r _(sh) +r _(hg) −r _(gs)

The subscripts stand as the Shuttle to Hubble, Hubble to a GPSsatellite, and A GPS satellite to the Shuttle respectively. The GPSsatellites position is known from the ephemeredes, and Shuttle'sposition is known from the receiver's navigation solution. If we define,Δ_(r)=X_(s)−X_(h), that is the relative position between the Hubble andShuttle and rotate into Hill's frame we arrive at:

ρ_(r)=|Δ_(r)|+|Δ_(r) −X _(g) |−|X _(g)|

Δ_(r) is our relative position vector to Hubble in the ECI frame, andpart of the state vector. Now taking partials with respect to Δ_(r) wearrive at:

$\frac{\partial\rho_{r}}{\partial\Delta_{r}} = {\frac{\Delta_{r}}{\Delta_{r}} + \frac{\Delta_{r} - X_{g}}{{\Delta_{r} - X_{g}}}}$$\frac{\partial\rho_{r}}{\partial{\overset{.}{\Delta}}_{r}} = 0$

This forms our measurement vector:

${H(x)} = \begin{bmatrix}\frac{\partial\rho_{r}}{\partial\Delta_{r\; 1x\; 3}} & 0_{1x\; 3}\end{bmatrix}$

The equation describing the reflected pseudorange rate measurement isgiven by:

h _({dot over (ρ)}r)(x)={dot over (ρ)}=(V _(s) −V _(h))·lôs_(sh)+(V _(h)−V _(g))·lôs_(hg)

Now, substituting in the line of sight vectors:

${\overset{.}{\rho}}_{r} = {{\left( {V_{s} - V_{h}} \right) \cdot \frac{X_{s} - X_{h}}{{X_{s} - X_{h}}}} + {\left( {V_{h} - V_{g}} \right) \cdot \frac{X_{h} - X_{g}}{{X_{h} - X_{g}}}}}$

Now substituting in Δ_(r)=X_(s)−X_(h), and Δ_(v)=V_(s)−V_(h):

$\overset{.}{\rho} = {{\Delta_{v} \cdot \frac{\Delta_{r}}{\Delta_{r}}} + {\left( {\Delta_{v} - V_{g}} \right) \cdot \frac{\Delta_{r} - X_{g}}{{\Delta_{r} - X_{g}}}}}$

Taking the partials with respect to the state Δr:

$\frac{\partial\overset{.}{\rho}}{\partial\Delta_{r}} = {{\Delta_{v}~ \cdot \left( {\frac{1}{\Delta_{r}} - \frac{\Delta_{r}^{2}}{{\Delta_{r}}^{3}}} \right)} + {\left( {\Delta_{v} - V_{g}} \right) \cdot \left( {\frac{1}{{\Delta_{r} - X_{g}}} + \frac{\Delta_{r} - X_{g}}{{{\Delta_{r} - X_{g}}}^{3}}} \right)}}$

Now, taking the partials with respect to the state Δv:

$\frac{\partial\overset{.}{\rho}}{\partial\Delta_{v}} = {\frac{\Delta_{r}}{\Delta_{r}} - \frac{\Delta_{r} - X_{g}}{{\Delta_{r} - X_{g}}}}$

Now to form the H vector:

${H_{\overset{.}{\rho}\; r}(x)} = \begin{bmatrix}\frac{\partial\overset{.}{\rho}}{\partial\Delta_{r}} & \frac{\partial\overset{.}{\rho}}{\partial\Delta_{v}}\end{bmatrix}$

For the nonlinear full dynamics implementation, measurements are givenby the measurement equations as a function of the states of the system,such that y=h(f). In this formulation, there are 4 differentmeasurements, pseudorange, differenced reflected pseudorange,pseudorange rate, and reflected pseudorange rate. In order to use thesemeasurements in our Kalman filter design, we must linearize themeasurements in order to form the canonical H matrix.The direct psuedorange measurement is given by the equation:

h _(ρ)(x)=ρ=r+c(t _(u) +δt+δt _(D))

neglecting the ionospheric, receiver, ephemeris and other errors(δ+δT_(d)), the measurement can expressed in terms of a range and bias.

${h_{\rho}(x)} = {\rho = {{r + b} = {\sqrt{\left( {x_{s} - x_{g}} \right)^{2} + \left( {y_{s} - y_{g}} \right)^{2} + \left( {z_{s} - z_{g}} \right)^{2}} + b}}}$

where b is the receivers time offset (in meters) from the GPSconstellation.

${H_{\rho}(x)} = \frac{\partial\rho}{\partial x}$

The partial with respect to the Shuttles position:

$\frac{\partial\rho}{\partial x_{s}} = \frac{x_{s} - x_{g}}{\left( {x_{s} - x_{g}} \right)^{2} + \left( {y_{s} - y_{g}} \right)^{2} + \left( {z_{s} - z_{g}} \right)^{2}}$${{let}\mspace{14mu} r_{sg}} = \sqrt{\left( {x_{s} - x_{g}} \right)^{2} + \left( {y_{s} - y_{g}} \right)^{2} + \left( {z_{s} - z_{g}} \right)^{2}}$

Now taking it with respect to each position vector:

$\begin{bmatrix}\frac{x_{s} - x_{g}}{r_{sg}} & \frac{y_{s} - y_{g}}{r_{sg}} & \frac{z_{s} - z_{g}}{r_{sg}}\end{bmatrix}$

The partial with respect to the Shuttle's velocity:

$\frac{\partial\rho}{\partial v_{s}} = 0_{1x\; 3}$

The partial with respect to receiver bias:

$\frac{\partial\rho}{\partial r_{b}} = 1$

The partial with respect to Hubble's position:

$\frac{\partial\rho}{\partial r_{h}} = 0_{1x\; 3}$

The partial with respect to Hubble's velocity:

$\frac{\partial\rho}{\partial v_{h}} = 0_{1x\; 3}$

Now forming a H vector for the direct measurement

${H_{\rho}(x)} = \begin{bmatrix}\frac{\partial\rho}{\partial r_{s\; 1x\; 3}} & 0_{1x\; 3} & 0_{1x\; 3} & 0_{1x\; 3} & 1\end{bmatrix}$

The differenced reflected pseudorange provides information about theabsolute shuttle and Hubble position states. Its measurement equation isgiven as follows:

h(x)=ρ_(r) =r _(sh) +r _(hy) −r _(sy)

The bias terms difference out, as well as many of the other common modeerrors greatly simplifying the calculation and the algebra. Thesubscripts sh, hg, and gs signify the range from shuttle to Hubble,Hubble to GPS transmitter, and shuttle to GPS transmitter. These areequivalently:

r _(sh)=√{square root over ((x _(s) ,x _(h))²+(y _(s) −y _(h))²+(z _(s)−z _(h))²)}{square root over ((x _(s) ,x _(h))²+(y _(s) −y _(h))²+(z_(s) −z _(h))²)}{square root over ((x _(s) ,x _(h))²+(y _(s) −y_(h))²+(z _(s) −z _(h))²)}

r _(hg)=√{square root over ((x _(h) −x _(g))²+(y _(h) −y _(g))²+(z _(h)−z _(g))²)}{square root over ((x _(h) −x _(g))²+(y _(h) −y _(g))²+(z_(h) −z _(g))²)}{square root over ((x _(h) −x _(g))²+(y _(h) −y_(g))²+(z _(h) −z _(g))²)}

r _(sg)=√{square root over ((x _(s) −x _(g))²+(y _(s) −y _(g))²+(z _(s)−z _(g))²)}{square root over ((x _(s) −x _(g))²+(y _(s) −y _(g))²+(z_(s) −z _(g))²)}{square root over ((x _(s) −x _(g))²+(y _(s) −y_(g))²+(z _(s) −z _(g))²)}

In order to form the H vector equation we must first linearize themeasurement equation. The partial with respect to the Shuttle'sposition:

$\frac{\partial\rho}{\partial x_{s}} = {\frac{x_{s} - x_{h}}{\sqrt{\begin{matrix}{\left( {x_{s} - x_{h}} \right)^{2} +} \\{\left( {y_{s} - y_{h}} \right)^{2} + \left( {z_{s} - z_{h}} \right)^{2}}\end{matrix}}} - \frac{x_{s} - x_{g}}{\sqrt{\begin{matrix}{\left( {x_{s} - x_{g}} \right)^{2} +} \\{\left( {y_{s} - y_{g}} \right)^{2} + \left( {z_{s} - z_{g}} \right)^{2}}\end{matrix}}}}$

substituting in the definitions of the range vectors we get:

$\frac{\partial\rho}{\partial x_{s}} = {\frac{x_{s} - x_{h}}{r_{sh}} - \frac{x_{s} - x_{g}}{r_{gs}}}$

Now extrapolating for each of the shuttle position components, we get a1×3 vector

$H_{x} = \begin{bmatrix}{\frac{x_{s} - x_{h}}{r_{sh}} - \frac{x_{s} - x_{g}}{r_{gs}}} & {\frac{x_{s} - x_{h}}{r_{sh}} - \frac{x_{s} - x_{g}}{r_{gs}}} & {\frac{x_{s} - x_{h}}{r_{sh}} - \frac{x_{s} - x_{g}}{r_{gs}}}\end{bmatrix}$

The partial with respect to the Shuttle's velocity:

$\frac{\partial\rho_{r}}{\partial{\overset{.}{x}}_{v}} = 0_{1x\; 3}$

The partial with respect to the Shuttle's bias:

$\frac{\partial\rho_{r}}{\partial{\overset{.}{x}}_{b}} = 0$

The partial with respect to Hubble's position:

$\frac{\partial\rho_{r}}{\partial x_{h}} = {{- \frac{x_{s} - x_{h}}{\sqrt{\begin{matrix}{\left( {x_{s} - x_{h}} \right)^{2} +} \\{\left( {y_{s} - y_{h}} \right)^{2} + \left( {z_{s} - z_{h}} \right)^{2}}\end{matrix}}}} + \frac{x_{h} - x_{g}}{\sqrt{\begin{matrix}{\left( {x_{h} - x_{g}} \right)^{2} +} \\{\left( {y_{h} - y_{g}} \right)^{2} + \left( {z_{h} - z_{g}} \right)^{2}}\end{matrix}}}}$

Now substituting in the range definitions, and taking partials of HST'sy and z components we arrive at the H vector:

$\frac{\partial\rho_{r}}{\partial x_{h}} = \begin{bmatrix}\begin{matrix}{{- \frac{x_{s} - x_{h}}{r_{sh}}} + \frac{x_{h} - x_{g}}{r_{hg}} -} \\{\frac{y_{s} - y_{h}}{r_{sh}} + \frac{y_{h} - y_{g}}{r_{hg}} -}\end{matrix} \\{\frac{z_{s} - z_{h}}{r_{sh}} + \frac{z_{h} - z_{g}}{r_{hg}}}\end{bmatrix}$

The partial with respect to Hubble's velocity:

$\frac{\partial{\overset{.}{r}}_{h}}{\partial\rho_{r}} = 0_{1x\; 3}$

Now forming a H vector for the direct measurement

${H_{\rho \; r}(x)} = \begin{bmatrix}\frac{\partial\rho}{\partial X_{s\; 1x\; 3}} & 0_{1x\; 3} & \frac{\partial\rho}{\partial X_{h\; 1x\; 3}} & 0_{1x\; 3} & 0\end{bmatrix}$

The equation that relates the Pseudorange rate measurement to the statesis given by:

${h(x)} = {\overset{.}{\rho} = {\left( {V_{u} - V_{g}} \right) \cdot \frac{X_{s} - X_{g}}{{X_{s} - X_{g}}}}}$

The measurement is simply the relative velocity along the line of sightvector between the GPS transmitter and the Shuttle.

${H_{\overset{.}{\rho}}(x)} = \frac{\partial\overset{.}{\rho}}{\partial X}$

This partials with respect to the scalar x position is given as:

$\frac{\partial\overset{.}{\rho}}{\partial x_{s}} = {\left( {V_{s} - V_{g}} \right){x\begin{bmatrix}{\frac{1}{\sqrt{\left( {x_{s} - x_{h}} \right)^{2} + \left( {y_{s} - y_{h}} \right)^{2} + \left( {z_{s} - z_{h}} \right)^{2}}} -} \\\frac{\left( {x_{s} - x_{g}} \right)^{2}}{\left( \sqrt{\left( {x_{s} - x_{h}} \right)^{2} + \left( {y_{s} - y_{h}} \right)^{2} + \left( {z_{s} - z_{h}} \right)^{2}} \right)^{3}}\end{bmatrix}}}$

Now substituting in the range formulas to clean up

$\frac{\partial\overset{.}{\rho}}{\partial x_{s}} = {\left( {V_{s} - V_{g}} \right)_{x}\left\lbrack {\frac{1}{r_{sh}} - \frac{\left( {x_{s} - x_{g}} \right)^{2}}{r_{sh}^{3}}} \right\rbrack}$

The same holds for the y and z vectors as well.

$\frac{\partial\overset{.}{\rho}}{\partial y_{s}} = {\left( {V_{s} - V_{g}} \right)_{y}\left\lbrack {\frac{1}{r_{sh}} - \frac{\left( {y_{s} - r_{g}} \right)^{2}}{r_{sh}^{3}}} \right\rbrack}$$\frac{\partial\overset{.}{\rho}}{\partial z_{s}} = {\left( {V_{s} - V_{g}} \right)_{z}\left\lbrack {\frac{1}{r_{sh}} - \frac{\left( {z_{s} - z_{g}} \right)^{2}}{r_{sh}^{3}}} \right\rbrack}$

Now taking the partials with respect to range, it simply leaves the lineof sight vectors.

$\frac{\partial\overset{.}{\rho}}{\partial V_{s}} = \frac{X_{s} - X_{h}}{{X_{s} - X_{h}}}$

breaking these into components and substituting in the range equation,we get:

$\frac{\partial\overset{.}{\rho}}{\partial V_{sx}} = \frac{x_{s} - x_{h}}{r_{s}h}$$\frac{\partial\overset{.}{\rho}}{\partial V_{sy}} = \frac{y_{s} - y_{h}}{r_{s}h}$$\frac{\partial\overset{.}{\rho}}{\partial V_{sz}} = \frac{z_{s} - z_{h}}{r_{s}h}$

Assuming a constant receiver bias:

$\frac{\partial\overset{.}{\rho}}{\partial X_{b}} = 0$

Because there is no information contained in the pseudorange rate aboutthe state vector of the Hubble, the partials with respect to Hubble'sstates are 0.

$\frac{\partial\overset{.}{\rho}}{\partial X_{h}} = 0_{1x\; 3}$$\frac{\partial\overset{.}{\rho}}{\partial V_{h}} = 0_{1x\; 3}$

Now constructing the H vector for the psuedorange rate measurement, weget:

${H_{\overset{.}{\rho}}(X)} = \begin{bmatrix}{\frac{\partial\overset{.}{\rho}}{\partial X_{s}}1x\; 3} & {\frac{\partial\overset{.}{\rho}}{\partial V_{s}}1x\; 3} & 0_{1x\; 3} & 0_{1x\; 3} & 0\end{bmatrix}$

The reflected pseudorange measurement is essentially the sum of theDoppler shift along the paths of the signal. This is given by:

h(X)={dot over (ρ)}=(V _(s) −V _(h))·lôs_(sh)+(V _(h) −V _(g))·lôs_(hg)

Now, expressing the line of sight vectors in terms of their components:

${h(X)} = {{\overset{.}{\rho}}_{rm} = {{\left( {V_{s} - V_{h}} \right) \cdot \frac{X_{s} - X_{h}}{{X_{s} - X_{h}}}} + {\left( {V_{h} - V_{g}} \right) \cdot \frac{X_{h} - X_{h}}{{X_{h} - X_{g}}}}}}$

From this, we need to take the partials with respect to the state vectorcomponents. For the shuttle position states, in the x direction:

$\frac{\partial\overset{.}{\rho}}{\partial x_{s}} = {\left( {V_{s} - V_{h}} \right)_{x}\left\lbrack {\frac{1}{{X_{s} - X_{h}}} - \frac{\left( {x_{s} - x_{h}} \right)^{2}}{{{X_{s} - X_{h}}}^{3}}} \right\rbrack}$

Substituting in the range vectors to simplify the equation we get

$\frac{\partial\overset{.}{\rho}}{\partial x_{s}} = {{\left( {V_{s} - V_{h}} \right)_{x}\left\lbrack {\frac{1}{r_{sh}} - \frac{\left( {x_{s} - x_{h}} \right)^{2}}{r_{hs}^{3}}} \right\rbrack}.}$

For the shuttle position states, in the y and z directions:

$\frac{\partial\overset{.}{\rho}}{\partial y_{s}} = {\left( {V_{s} - V_{h}} \right)_{y}\left\lbrack {\frac{1}{r_{sh}} - \frac{\left( {y_{s} - y_{h}} \right)^{2}}{r_{hs}^{3}}} \right\rbrack}$$\frac{\partial\overset{.}{\rho}}{\partial z_{s}} = {\left( {V_{s} - V_{h}} \right)_{z}\left\lbrack {\frac{1}{r_{sh}} - \frac{\left( {z_{s} - z_{h}} \right)^{2}}{r_{hs}^{3}}} \right\rbrack}$

Now taking the partial with respect to the Shuttle's velocity:

$\frac{\partial\overset{.}{\rho}}{\partial V_{s}} = \frac{X_{s} - X_{h}}{{X_{s} - X_{h}}}$

This simply leaves the line of sight vector to the Hubble; expandingthis into the x,y, and z components:

$\mspace{79mu} {\frac{\partial\overset{.}{\rho}}{\partial V_{sx}} = \frac{x_{s} - x_{h}}{r_{sh}}}$$\mspace{79mu} {\frac{\partial\overset{.}{p}}{\partial V_{sy}} = \frac{y_{s} - y_{h}}{r_{sh}}}$$\mspace{79mu} {\frac{\partial\overset{.}{p}}{\partial V_{sz}} = \frac{z_{s} - z_{h}}{r_{sh}}}$$\frac{\partial\overset{.}{\rho}}{\partial X_{hx}} = {{{- \left( {V_{s} - V_{h}} \right)_{x}} \cdot \left\lbrack {\frac{1}{{X_{s} - X_{h}}} - \frac{\left( {x_{s} - x_{h}} \right)^{2}}{{{X_{s} - X_{h}}}^{3}}} \right\rbrack} + {\left( {V_{h} - V_{g}} \right) \cdot \left\lbrack {\frac{1}{{X_{h} - X_{g}}} - \frac{\left( {x_{h} - x_{g}} \right)^{2}}{{{X_{h} - X_{g}}}^{3}}} \right\rbrack}}$

substituting in the range vectors to clean up the equation, andextrapolating to the y,z vectors we arrive at:

$\frac{\partial\overset{.}{\rho}}{\partial X_{hx}} = {{{- \left( {V_{s} - V_{h}} \right)_{x}} \cdot \left\lbrack {\frac{1}{r_{sh}} - \frac{\left( {x_{s} - x_{h}} \right)^{2}}{r_{sh}^{3}}} \right\rbrack} + {\left( {V_{h} - V_{g}} \right)_{x} \cdot \left\lbrack {\frac{1}{r_{hg}} - \frac{\left( {x_{h} - x_{g}} \right)^{2}}{r_{hg}^{3}}} \right\rbrack}}$$\frac{\partial\overset{.}{\rho}}{\partial X_{hy}} = {{{- \left( {V_{s} - V_{h}} \right)_{y}} \cdot \left\lbrack {\frac{1}{r_{sh}} - \frac{\left( {y_{s} - y_{h}} \right)^{2}}{r_{sh}^{3}}} \right\rbrack} + {\left( {V_{h} - V_{g}} \right)_{y} \cdot \left\lbrack {\frac{1}{r_{hg}} - \frac{\left( {y_{h} - y_{g}} \right)^{2}}{r_{hg}^{3}}} \right\rbrack}}$$\frac{\partial\overset{.}{\rho}}{\partial X_{hy}} = {{{- \left( {V_{s} - V_{h}} \right)_{z}} \cdot \left\lbrack {\frac{1}{r_{sh}} - \frac{\left( {z_{s} - z_{h}} \right)^{2}}{r_{sh}^{3}}} \right\rbrack} + {\left( {V_{h} - V_{g}} \right)_{z} \cdot \left\lbrack {\frac{1}{r_{hg}} - \frac{\left( {z_{h} - z_{g}} \right)^{2}}{r_{hg}^{3}}} \right\rbrack}}$

Now taking the partial with respect to the Hubble's velocity, we simplyget a relation involving the line of sight vectors.

$\frac{\partial\overset{.}{p}}{\partial V_{h}} = {{- \frac{X_{s} - X_{h}}{{X_{s} - X_{h}}}} + \frac{X_{h} - X_{g}}{{X_{h} - X_{g}}}}$

expanding it into the x,y,z components and substituting in the rangedefinitions we get:

$\frac{\partial\overset{.}{\rho}}{\partial V_{hx}} = {{- \frac{x_{s} - x_{h}}{r_{sh}}} + \frac{x_{h} - x_{g}}{r_{hg}}}$$\frac{\partial\overset{.}{\rho}}{\partial V_{hy}} = {{- \frac{y_{s} - y_{h}}{r_{sh}}} + \frac{y_{h} - y_{g}}{r_{hg}}}$$\frac{\partial\overset{.}{\rho}}{\partial V_{hz}} = {{- \frac{z_{s} - z_{h}}{r_{sh}}} + \frac{z_{h} - z_{g}}{r_{hg}}}$

Since the reflected pseudorange rate does not involve the bias, it iszero.

$\frac{\partial\overset{.}{\rho}}{\partial X_{b}} = 0$

Forming the linearized measurement vector we get:

${H_{\overset{.}{\rho \; r}}(X)} = \left\lbrack {\frac{\partial\overset{.}{\rho \; r}}{\partial X_{s\; 1 \times 3}}\mspace{14mu} \frac{\partial\overset{.}{\rho \; r}}{\partial V_{s\; 1 \times 3}}\mspace{14mu} \frac{\partial\overset{.}{\rho \; r}}{\partial X_{h\; 1 \times 3}}\mspace{14mu} \frac{\partial\overset{.}{\rho \; r}}{\partial V_{h\; 1 \times 3}}\mspace{14mu} 0} \right\rbrack$

The measurement matrix corrects the Kalman Filter's dynamic estimateswith the current measurement estimates. The measurement matrix isstacked up in the same order as the measurements.The measurement for the relative filter is given by:

H=[H _(ρr)(X),H _({dot over (ρ)}r)(X)]^(T)

The matrix for the Absolute dynamic filter is given as:

H=[H _(ρ)(X),H _({dot over (ρ)})(X),H _({dot over (ρ)}r)(X),H_({dot over (ρ)}r)(X)]^(T)

The process noise matrix, Q describes the covariance error associatedwith the model. It generally describes the model errors and un-modeledforce inputs. In our case, since Q is difficult to obtain, it is treatedas a tuning matrix, to adjust the performance of the filter.

$Q = \begin{bmatrix}{\sigma_{p}^{2}I_{3 \times 3}} & 0_{3 \times 3} \\0_{3 \times 3} & {\sigma_{v}^{2}I_{3 \times 3}}\end{bmatrix}$

σ_(p) is the standard deviation in Hill's relative position model, σ_(v)is the standard deviation in Hill's relative velocity model.

$Q = \begin{bmatrix}{\sigma_{ps}^{2}I_{3 \times 3}} & 0 & 0 & 0 & 0 \\0 & {\sigma_{vs}^{2}I_{3 \times 3}} & 0 & 0 & 0 \\0 & 0 & {\sigma_{ph}^{2}I_{3 \times 3}} & 0 & 0 \\0 & 0 & 0 & {\sigma_{vh}^{2}I_{3 \times 3}} & 0 \\0 & 0 & 0 & 0 & \sigma_{b}^{2}\end{bmatrix}$

σ_(ps) is the standard deviation in the Shuttle's position model, σ_(vs)is the standard deviation in the Shuttle's velocity model. σ_(ph) andσ_(vh) represent the standard deviation in Hubble's position andvelocity models respectively.

The R matrix represents the variance in the measurements. Ideally thevariance is uncorrelated white gaussian noise. Since we can estimate thevariance of the measurements with relative accuracy, this matrix isessentially fixed.

$R = \begin{bmatrix}{\sigma_{\rho \; r}^{2}I_{m \times m}} & 0_{m \times m} \\0_{m \times m} & {\sigma_{\overset{.}{\rho \; r}}^{2}I_{m \times m}}\end{bmatrix}$

σ_(Dρ) and σ_(D{dot over (ρ)}) are the differenced psuedorange andreflected pseudorange rate measurements.

$R = \begin{bmatrix}{\sigma_{\rho}^{2}I_{n \times n}} & 0 & 0 & 0 \\0 & {\sigma_{\overset{.}{\rho}}^{2}I_{n \times n}} & 0 & 0 \\0 & 0 & {\sigma_{\rho \; r}^{2}I_{m \times m}} & 0 \\0 & 0 & 0 & {\sigma_{\overset{.}{\rho \; r}}^{2}I_{n \times n}}\end{bmatrix}$

σ_(ρ) is the standard deviation in the pseudorange measurement,σ_({dot over (ρ)}) is the standard deviation in the pseudorange ratemeasurement. σ_(Dρ) and σ_(D{dot over (ρ)}) are the differencedpsuedorange and reflected pseudorange rate measurements.

Implementation of a preferred embodiment will now be described.

The use of reflected GPS signals as a navigational tool may beparticularly useful when navigating toward an orbiting target such asthe Hubble Space Telescope (HST) in order to perform repairs, replacecritical parts or otherwise service an orbiting device. The ability topassively and autonomously navigate to perform servicing missions couldlead to unmanned service missions as well as reduced weight and powerconsumption. The following description relates to navigating the spaceshuttle to the HST.

The space shuttle will incorporate a GPS receiver such as that disclosedand describe in US 2006-0082496 to Winternitz et al. which isincorporated herein in its entirely. The '496 receiver is particularlysuited for implementing the method of the present invention usingreflected GPS signals by having the ability to acquire and track weakGPS signals. The navigation system will include two antennas; a RightHanded and Left Handed Circular Polarized (LHCP, RHCP) antennas. The GPSreceiver will be powered on as it approaches the HST within 1000 m fromHubble. The Shuttle's position can be estimated with conventionalmethods of navigation using direct GPS signals and a dynamic orbitmodel. The receiver onboard the Space Shuttle will be acquiring thedirect signals as it approaches the HST which will be received by theRHCP GPS antenna. A convenient feature about a RHCP signal is that, whenreflected off of a surface such as the HST, it reverses polarization andbecomes left hand circularly polarized. This gives it a nearly perfectrejection for a reflected signal. Therefore, the LHCP Antenna is wellsuited to pick up the reflected signals.

Both of the antennae will receive the raw RF stream mixed down to IF,35.42 MHz, and then sampled at 8.192 MS/s. This data will be processedas previously described to compare the direct and reflected signals toestimate the position of the Shuttle relative to the HST.

The use of ever increasing speeds of computers, microprocessor andFPGAs, complex signal processing algorithms will be used to apply theaforementioned filtering techniques and signal acquisition methods todetermine the structure and composition of the reflected signals.

The GPS receiver used to record and process the data is the NavigatorGPS. The GPS receiver is running at 65.536 MHz using a General DynamicsColdfire processor and has three application specific Actel AX2000FPGA's. The enhanced capabilities of this receiver are aided by theapplication specific FPGAs that implement an algorithm developed toacquire weak signals in real time as disclosed '496 to Winternitz et al.This algorithm is an FFT based acquisition algorithm that operates inthe frequency domain. It is able to acquire and track signals down to 25dB-Hz. This is approximately 10 dB more sensitive than current off theshelf GPS receivers. This will enable the Navigator to acquire theweaker reflected signals off of the HST.

The Navigator receiver will be configured to track 24 channels, 12 ofthem connected to the LHCP antenna and the remaining 12 connected to theRHCP antenna. While tracking a direct signal, the reflected channelswill be directed to search in a time delayed space with approximatelythe same Doppler frequency. This will tremendously aid acquisition ofthe weaker reflected signals.

The aforementioned extended Kalman filter will run on the navigationprocessor taking in the RHCP and LHCP measurements. Thus the Navigationsystem will be able to estimate the position of the Shuttle in orbit andthe relative position of the HST by comparing the directed and reflectedsignals of a plurality of GPS satellites as previously discussed. Usinga receiver to receive both direct and reflected GPS signals andcomparing/processing the embedded ranging data (time delay etc.),coupled with the dynamic model in the extended Kalman filter, thenavigation system of the present invention can passively andautonomously track a relative position of an orbiting spacecraft(Shuttle) to an orbiting target (HST). The ability to passively andautonomously passively navigate a space craft toward a target not onlyreduces weight and power consumption, but also paves the way forpotential precise passive navigation of space craft for unmanned servicemissions.

A simulation and implementation of the method of the present inventionwill now be discussed in regards to a Shuttle mission approaching theHubble Telescope.

The Hubble Rendezvous truth-scenarios are chosen such that the initialconditions will result in a near collision of the Shuttle and Hubble.The truth model is generated with STK using a high fidelity orbitpropagator.

The Space Shuttle may approach Hubble from the aft end with its baydoors pointed towards Hubble. This allows the Shuttle to dock cleanlywith Hubble.

During the servicing mission, the left and right handed GPS Antennaewill be placed on the MULE inside the bay of the shuttle lookingoutwards.

The model focuses on simulating the uncontrolled motion of the twobodies. In order to represent the “docking”, initial conditions werechosen such that the two space craft come within close proximity to eachother.

The model was created using STK, given a set of initial conditions, andpropagated over an orbit period. The STK model chosen was a 21-ordergravity model with atmospheric drag on both the Shuttle and Hubble, thecoefficients of drag chosen as to represent the relative mass and arearatios.

Units Shuttle ECI Position (−4673.7, −4525.4, 2443.5)^(τ) Km ECIVelocity  (5.604, −4.472, 2.438)^(τ) Km/s Mass 109,000 kg Area 58.4 m²C_(d) 1.2 Hubble ECI Position (−4673.9, −4525.9, 2443.7)^(τ) Km ECIVelocity (−5.604, −4.471, 2.437)^(τ) Km/s Mass 11,110 kg Area 15.4 m²C_(d) 0.8The area for the Shuttle was approximated using the circular crosssectional area of the Shuttle. The same was done for the Hubble.

The visibility of each GPS satellite is determined by the geometry ofthe GPS constellation with respect to the location of Hubble and theShuttle. The satellite visibility is essential to determining how manymeasurements both direct and reflected are available. This number ofmeasurements available influences the absolute and relative positionaccuracies.

The GPS constellation is constructed and propagated using the keplerianelements given by an almanac file. The almanac file is generated anduploaded to each GPS satellite. It represents a limited set of thekeplerian orbits, in order to provide a rough position estimate of theall of GPS satellites. These orbits are then propagated in time,providing a good estimate of the satellites orbits. This is purely usedto generate a moving rough estimate of the GPS constellation for thesimulation.

In practice the GPS constellation is closely monitored, and eachsatellite has a new ephemeris updated every 3 hours. In addition toupdating their GPS Almanac files every 12 hours.

The Almanac parameters are given by:

Parameter Symbol Scaling Time of Almanac toa secs Eccentric Anamoly Enone Inclination i0 semicircles + 0.94247 Rate of Right Ascention {dotover (Ω)} semicircles Semi-major axis A½ none longitude of perigee λsemicircles Argument of Perigee ω semicircles Mean Anamoly M semicircles

Since we are only interested in observing those GPS Satellites with areflected pair, and since the relative distance is so small, we canassume that any GPS satellite seen by the space Shuttle can also be seenby the Hubble.

The visible signals are done by masking out those GPS satellites whosemain lobes cross through the earth, or are not in view of the shuttle.This is done by calculating the angle between the line of sight vector,and position vector of the shuttle as seen in FIG. 4. If β is greaterthan 13.9° and less than 21.3° then it is within the main lobe and justover the lim of the earth. If β is less than 13.9°, then we need tocheck to make sure the signal does not pass through the earth. Thus, thedot product of the two vectors should be negative.

The bi-static radar equation described in is given as:

$\left( {R_{T}R_{R}} \right) = \left( \frac{P_{T}G_{T}G_{R}\lambda^{2}\sigma_{B}F_{T}^{2}F_{R}^{2}}{\left( {4\pi} \right)^{2}{KT}_{S}{B_{N}\left( {S\text{/}N} \right)}_{\min}L_{T}L_{R}} \right)^{1/2}$R_(T) Transmitter to target range R_(R) Target to receiver range P_(T)Transmitter power G_(T) Transmit antenna power gain G_(R) Receiveantenna power gain σ_(b) Bistatic radar cross section F_(T) Patternpropagation factor for transmit to target F_(R) Pattern propagationfactor for receiver to target K Boltzman's constant T_(n) Receiver noisetemperature D_(n) Noise bandwidth of the receivers pre-detection filter(S/N) signal-to-noise ratio required for detection L_(T) Transmit systemlosses L_(R) Receive system losses

This equation is used to determine the maximum range to the target.However, rearranging terms, it is possible to pose this in the form ofSNR received.

There are several universal constants from this particularimplementation. λ is the wavelength of the GPS L1 signal. The GPS L1signal is transmitted at 1575.42 MHz, giving it a wavelength of 0.1903m. Boltzmann's constant is given by $1.3806503×10⁻²³ J/K. F_(T) andF_(R) are the pattern propagation factors, these factors describe theloss due to transmission in an unclean environment. F_(T) and F_(R) canbe estimated at around 2 dB total for the ionospheric effects.

The bi-static radar cross section is calculated in much the same way asthe standard radar cross section. However, the transmit angle needs tobe taken into account, as well as properties limited only to thebi-static scenario.

Since the Hubble can be approximated as a cylinder and the approachprofile is towards the aft-end. The Hubble can approximated by thereflective properties of a circular flat plate combined with acylindrical surface normal to the flat plate. In physical optics, theequation describing the RCS for a metallic circular disk is:

$\sigma = {16\; \pi {{\frac{A\; \cos \; \theta}{\lambda} \cdot \frac{J_{1}\left( {{kd}\; \sin \; \theta} \right)}{{kd}\; \sin \; \theta}}}^{2}}$

Where A is given as the physical area of the disk, d is the diameter,and J₁(x) is the Bessel function of the first kind of order one. λ isthe wavelength.

The Hubble space telescope has a diameter of roughly 4.2 m. Thus givingan area of roughly 13.85 m².

For normal incidence, this equation simplifies to:

$\sigma = \frac{4\; \pi \; A^{2}}{\lambda^{2}}$

This would represent the monostatic cross sectional area of the Hubblespace telescope with the observer looking perpendicular to the aft end.Thus, the maximum possible cross-sectional area, looking at the aft endis 36.81 dB-sqm. It is interesting to note that the GPS signal will getgain from the reflections. If the observer was oriented perpendicular toHubble, it would appear as a cylinder to the viewer. In this case, theRCS is given by:

$\sigma = {{kal}^{2}{\frac{\sin \left( {{kl}\; {\sin ({theta})}} \right)}{{kl}\; {\sin ({theta})}}}^{2}}$

Where α=2.1 m is the radius, l=14.2 m is the length of the Hubble, k isthe wave number, and λ is the wavelength. This provides a predictableradar cross section. With the assumption since at the extremes of bothof these figures, the RCS drops off significantly, and effectively takethe maximum RCS. Of the two plots, the side of the cylinder offset by90° to form a combined RCS model for the cylinder over 180°. Since acylinder is symmetrical, this is effectively models the full cylinderover 360°.

In order to account for the bi-static nature of the problem, it isnecessary to define a new angle called the bi-static angle. This angleis defined as the difference between the receive and transmit antennae,β=θ_(r)−θ_(t), the bisection of this angle is used to approximate theBi-static RCS, σ_(b). The bi-static RCS can be grouped into 3 sections:the pseudo-monostatic RCS region, the forward scatter RCS region, andthe bi-static RCS region. In the pseudo-monostatic region, the Crispinand Seigal monostatic-bi-static equivalence theorem predicts that for asufficiently small angle, the RCS is equivalent to the monostatic RCSmeasured the bisector of the bistatic angle. This theorem applies to asufficiently smooth surface, such as a sphere, a cylinder, or cone. Someempirical data for a cylinder, performed at 35 GHz, thepseudo-monostatic region extends to roughly 20°, while the bistatic iswithin 20°≦β≦140° and forward scatter at β>140 °.

The bi-static region begins where the pseudo-monostatic region ends, itis where the theorem fails to predict an accurate RCS. Unfortunately,there is no empirical formula for this equation. Experimental data showsa downward trend, from −2 dB to −15 dB, in this case. It should sufficeto use a linear approximation (in dB).

The third bi-static region is the forward scattered region, this occurswhen the bi-static angle approaches $180̂o$. The forward scatter can beapproximated by treating the shadow area as a uniformly illuminatedantenna aperture. The roll off can be approximated when {dot over (Π)}-βis substituted for the angle off aperture normal. This function closelymatches J₀(x)/x where J₀(x) is the Bessel function of order 0. In oursituation, many of the GPS satellites will lie in this area, due tobeing below the constellation. This is to our benefit, allowing us topick up signals, being at such an odd orientation. Luckily due to thechange in polarization of the reflected signal, this will provide uswith isolation from the direct signal from the satellite.

From these bi-static regions we can now form a full approximation forthe bi-static radar coefficient. The plot below represents σ_(b) withthe flat face of the cylinder oriented towards the receiver. This shouldprovide a reasonable fit for a perfectly reflecting disk.

The antennae chosen for this mission are two hemispherical antennas. Theantenna gain is given as 3 dBic. Built into each antenna is the lownoise amplifier. The LNA has a gain of 26 dB with a noise figure of 2.8.This can be used to calculate the system noise temperature. Thisequation is given as:

T _(s)=10 log(T _(A)(NF−1))

Where, T_(A) is the temperature at the antenna. Since the LNA is at theantenna, the line losses are included in the analysis. From the specsheet, the NF of the antenna is given as 2.8 dB. A conservativeapproximation for this noise figure is to use the temperature at theantenna, T_(A), as room temperature or 290 K. This approximation holdsbecause the antennae are lying on the MULE. The MULE will be heated inorder to keep the electronics from freezing. From this equation we findthat the receiver noise temperature is 24.4 dB-K.

The L_(R) and L_(T) and receiver and transmit system losses respectivelywill be taken into account as a total system loss. A generalapproximation of loss in a communications system is about 2 dB. Usuallythis value is measured quantitatively, but, without a complete system totest, it is difficult to arrive at a true figure. The noise bandwidth,B_(n) is generally taken as 1/T, where T is the pre-integration period.For Navigator, it has an adjustable preintegration period, however thenominal case is 1 ms. Other losses occur such as Polarization mismatchcontributing 3 dB and in the case of a reflection, there is anefficiency which is dependent on the reflectivity of the object. In thecase of the Hubble, which is coated in mylar, it will have a highreflectivity coefficient, this is typically 90-98%, depending on thequality.

The range calculations R_(R) is simply the geometric distance betweenthe Shuttle and the Hubble, while R_(T) is simply the geometric distancebetween each GPS Satellite and the Hubble. Factoring the free space lossdue to range the equations:

${R_{T}({dB})} = {20\; {\log \left( \frac{4\pi \; R_{T}}{\lambda} \right)}}$${R_{R}({dB})} = {20\; {\log \left( \frac{4\pi \; R_{R}}{\lambda} \right)}}$

At each time-step, for each visible GPS satellite, each range value iscalculated. With the assumption that the Shuttle is oriented with itsbay doors towards the aft bulkhead of Hubble, it is possible to estimatethe angle β and which surface of the cylinder the reflection isstriking. The truth positions were generated from a file that is thesimple 2-body relative motion, coupled with some thruster firing for arealistic closing profile. This produces an expected satellitevisibility given the aforementioned conditions. It should be noted thatthe model does not take into account the reflections from the solarpanels and assumes an orientation with the antennae directly facing theaft bulkhead.

The direct signal is generated using the truth positions of the Shuttleand Hubble. A receiver bias is added on top of the signal, as well asthe ionospheric error given by the equation

$\delta_{iono} = \frac{82.1\; T\; E\; C}{{F_{C}^{2}\sqrt{{\sin^{2}({EL})} + 0.76}} + {\sin^{2}({EL})}}$

TEC is the total electron content, given by 2e17. F_(c) is the carrierfrequency 1.57542 GHz, and EL is the elevation of the GPS satellite overthe horizon.

The direct signal pseudorange rate is generated by the rate of change ofthe line of sight vectors. In practice, most GPS receivers rely only onthe pseudorange measurements or use carrier phase to generate thevelocities. The pseudorange rate measurements are generally noisy andactually degrade the performance of the filter.

A variance of 5 meters was chosen for the direct pseudo-ranges, and 25cm/s for the pseudo-range rates. It is composed of the clock drift,ephemeris errors, phase noise, and other unmodeled sources. Theionospheric bias was neglected in this simulation, as it would have beendifferenced out in the measurement generation. The pseudo-range variancemay be increased by a nominal amount, 0.25 m, to reflect the ionosphericmodel discrepancies.

The differenced reflected measurements were constructed from the totalpath lengths between the GPS Satellite, Hubble, and the Shuttle. Thenoise added is from the loss in signal quality due to the reflectiveproperties and small ionospheric error. Since the measurement is adifferenced reflected signal with that of the direct path, the bias andthe ionospheric terms become small and are negligible. In addition tothe ionospheric error practically canceling out, the front end bias,clock bias, thermal jitter and other common mode errors also disappear.

A variance of 25 meters was chosen for the reflected pseudo-ranges, and50 cm/s for pseudo-range rates. Although the true variance is unknown,it is extremely difficult to model on the earth. Future tests willinclude an anechoic chamber testing or data reduction of the rendezvousexperiment.

While the method and apparatus of the present invention has yet to betested oil a service mission, simulations have proven the viability andutility of the use of comparing reflected and directly received GPSsignals to passively navigate spacecraft to an orbiting target. Thepresent invention is currently scheduled to fly in August 2008 on HubbleServicing Mission 4. A further discussion of such simulations can befound in Ian Cohen's thesis entitled “Relative Navigation for HubbleMission Using Reflected GPS signals”, published by the University ofMaryland in May 2007 and is hereby incorporated herein by reference inits entirety.

One additional use of the method of using reflected GPS signals as apassive navigation system is to compliment existing navigation tools.This passive system can not only provide an extra measurement tooloffering additional redundancy, it could also be used to navigate aspace craft to within sufficient proximity of an orbiting targetreducing the amount of time active systems such as LIDAR and RADAR arepowered. Because this method and system of the present invention may beemployed by simply adding the second LHCP antenna and loading thealgorithmic complexity onto existing systems, redundancy and reducedpower consumption can be achieved with little modifications to existingequipment.

The foregoing description of the specific embodiments will so fullyreveal the general nature of the invention that others can, by applyingcurrent knowledge, readily modify and/or adapt for various applicationssuch specific embodiments without departing from the generic concept,and, therefore, such adaptations and modifications should and areintended to be comprehended within the meaning and range of equivalentsof the disclosed embodiments. It is to be understood that thephraseology or terminology employed herein is for the purpose ofdescription and not of limitation. Therefore, while the embodiments ofthe invention have been described in terms of preferred embodiments,those skilled in the art will recognize that the embodiments of theinvention can be practiced with modification within the spirit and scopeof the appended claims.

1. A method of passively navigating an orbiting spacecraft towards anorbiting target, said method comprising the steps of: receiving aplurality of direct signals, each of said plurality of direct signalstransmitted from a corresponding one of a plurality of celestial GPSearth orbiting satellites transmitting ranging data signals; receiving areflected signal from said corresponding one of said plurality of GPSsatellites reflected off said orbiting target; comparing said directsignals and said reflected signals; and determining a position of saidorbiting spacecraft relative to said orbiting target based on said stepof comparing said direct signals and said reflected signals.
 2. Themethod of claims 1, wherein said step of comparing said direct andreflected signals includes the step of determining a time delay betweenreceiving said direct signal and receiving said reflected signal from aparticular GPS satellite.
 3. The method according to claim 1 whereinsaid steps of receiving said direct and reflected and said step ofcomparing said direct signals and said reflected signals areautonomously and passively performed on said orbiting space craft. 4.The method of claim 1, step of determining a relative position of saidorbiting spacecraft includes tracking the position of the spacecraft andupdating the position of the spacecraft as a function of the trackedposition and said ranging data.
 5. The method of claim 4, wherein stepof determining said relative position of said orbiting spacecraftincludes filtering said ranging data in said reflected signals and saiddirect signals.
 6. The method of claim 5, wherein step of determiningsaid relative position of said space craft includes employing anextended Kalman filter to compute the position of the orbitingspacecraft relative to said orbiting target.
 7. The method according toclaim 1, wherein said step of receiving said direct signals is performedwith a left handed circular polarized antenna and said step of receivingsaid reflected signal is performed with a right handed circularpolarized antenna.
 8. The method according to claim 1, wherein said stepof comparing said direct and reflected signals includes processing saidsignals into pseudorange measurements.
 9. The method of claim 8, whereinsaid step of processing includes processing said direct signal into apseudorange measurement representing a distance directly between saidorbiting space craft and a corresponding satellite and processing saidreflected signals into a pseudorange measurement representative of a sumdistance between said orbiting space craft to said orbiting target andsaid orbiting target to said corresponding satellite.
 10. The method ofclaim 3, further comprising a step of implementing a extended Kalmanfilter to process said ranging data to generate and update estimatedposition data thereby autonomously generate real time estimates of saidrelative position of said orbiting space craft relative to said orbitingtarget.
 11. A navigation system to passively navigate an orbitingspacecraft to an orbiting target, said navigation system comprising: areceiver platform adapted to receive a plurality of signals carryingranging data transmitted by at least one of a constellation ofsatellites; said receiving platform including; a first antenna adaptedto receive at least one of said plurality of signals direct from acorresponding satellite; a second antenna adapted to receive said atleast one of said plurality of signals from said corresponding satellitereflected off said orbiting target; a processor adapted to process saidranging data received by said first and second antennas to therebyautonomously and passively locate a relative position of said orbitingspacecraft relative to said orbiting target.
 12. The navigation systemaccording to claim 1 further comprising; an extended Kalman filteradapted to process said ranging data received by said receiver platformand generate and update estimated position data to autonomously generateestimates of said relative position of said orbiting space craftrelative to said orbiting target.
 13. The navigation system according toclaim 1, wherein said processor is adapted to process said at least oneof said plurality of signals direct from a corresponding satellite intoa pseudorange measurement representing a distance directly between saidorbiting space craft and said corresponding satellite and processingsaid at least one of said plurality of signals reflected of saidorbiting target into a pseudorange measurement representative of a sumdistance between said orbiting space craft to said orbiting target andsaid orbiting target to said corresponding satellite.
 14. The passivenavigation system according to claim 1, wherein said first antenna is aleft handed circular polarized antenna and said second antenna is aright handed circular polarized antenna.
 15. The passive navigationsystem according to claim 4, wherein each of said first and secondantennas is a hemispherical antenna each incorporating a low noiseamplifier.
 16. The passive Navigation system according to claim 13,wherein said receiving platform receives signals transmitted from atleast four satellites and said processor implements said extended Kalmanfilter to measure a relative position of said orbiting spacecraftrelative to said orbiting target.
 17. A method of passively navigating amoving body a moving target, said method comprising the steps of:receiving a plurality of direct signals, each of said plurality ofdirect signals transmitted from a corresponding one of a plurality ofcelestial GPS earth orbiting satellites transmitting ranging datasignals; receiving a reflected signal from said corresponding one ofsaid plurality of GPS satellites reflected off said moving target;implementing an extended Kalman filter to process said ranging data togenerate and update estimated position data thereby autonomouslygenerate real time estimates of said relative position of said movingspace craft relative to said moving target; comparing said directsignals and said reflected signals; and determining a position of saidmoving body relative to said moving target based on said step ofcomparing said direct signals and said reflected signals.